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Abstract 

We compute properties of the interface of the 3-dimensional Ising model 
for a wide range of temperatures, covering the whole region from the low tem- 
perature domain through the roughening transition to the bulk critical point. 
The interface tension a is obtained by integrating the surface energy density 
over the inverse temperature (3. We use lattices of size L x L x T, with L 
up to 64, and T up to 27. The simulations with antiperiodic boundary con- 
ditions in T-direction are done with the Hasenbusch-Meyer interface cluster 
algorithm that turns out to be very efficient. We demonstrate that in the 
rough phase the large distance behavior of the interface is well described by a 
massless Gaussian dynamics. The surface stiffness coefficient k is determined. 
We also attempt to determine the correlation length £ and study universal 
quantities like £ 2 <r and £ 2 k. Results for the interfacial width on lattices up to 
512 x 512 x 27 are also presented. 



1 Introduction 



There has been continuous interest in the properties of interfaces separating coex- 
isting phases. A prominent role is played by the 3- dimensional Ising model that is 
believed to share a universality class with binary systems in nature. The dominating 
method for quantitative studies of the interface of the 3-dimensional Ising model is 
the Monte Carlo method [|IJ, see e.g. 0, [| |, || || [7], and references cited therein. 

The Ising interface undergoes a roughening transition at an inverse temperature 
(3 T that is nearly twice as large as the bulk transition coupling (3 C (= 0.221652(3) 
U). The most precise estimate for the roughening coupling is /3 r = 0.4074(3) |T0| . 



For (3 C < (3 < (3 r the interface is rough. It is believed that its infrared properties 
can be described by a massless Gaussian dynamics. This is the basic assumption of 
the theory of capillary waves that is widely used to describe long distance properties 
of rough interfaces. A Gaussian behavior in the infrared is also what is expected 



from a Kosterlitz-Thouless model in the massless phase |TTJ. Indeed, it is believed 
that the roughening transition of the 3-dimensional Ising model is of the Kosterlitz- 
Thouless type. This belief is strongly substantiated by the renormalization group 



analysis in |10 



At and above the bulk transition temperature the system properties become 
independent of the boundary conditions in the infinite volume limit. The interface 
tension (free energy per unit area of the interface) vanishes like a ~ oot^, with 
t — \{(3 — (3 C )/(3 C \. If Widom scaling holds, the exponent /i should be twice the 
exponent v that determines the critical behavior of the correlation length. The 
most precise estimates for v are in the range to 0.624.. .0.630 |, [12|. The amplitude 
cr is of particular interest because it enters certain universal amplitudes that can 
also be measured in real life systems. 

In this paper, we report on a numerical study of properties of the Ising interface 
over a wide range of temperatures: from the low temperature regime through the 
rough domain up to the bulk transition region. The focus is mainly on the interface 
free energy, the interface tension, and on the long distance properties in the rough 
phase. 

The interface free energy is determined by integrating the surface energy over f3. 
(To the best of our knowledge, this method to obtain surface free energies was first 
used by Burkner and Stauffer |13fl ). We start the integration both at the high and 
at low temperatures and compare the results. 

Our method allows to include interfaces with extension up to 64 x 64 in the 
analysis. Close to criticality we also use a finite step method that allows to directly 
obtain the change of the surface free energy over a small interval Aj3. 

From the interface free energies we get estimates for the interface tension a and 
make fits with the critical law cited above. We also determine the correlation length 
£ in order to study the quantity £ 2 <r where the factors t~ 2u and t** cancel each other. 

In the rough phase, we study the long distance behavior of the interface by 
measuring block spin correlation functions of suitably defined interface "height vari- 
ables". An effective coupling /3 e ff (well known in Kosterlitz-Thouless theory) is 
obtained that parameterizes the asymptotic Gaussian dynamics. /3 e s is related to 
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the surface stiffness coefficient k that enters the surface Hamiltonian of the capillary 
wave model. 

We also study the surface width. In the rough phase, the width is expected to 
grow logarithmically with the surface extension, with a coefficient that is propor- 
tional to ft e g. For a simulation at ft — /3 c /0.8 we verify this behavior with good 
accuracy. 

The paper is organized as follows: In section |2|, we introduce the model and our 
conventions. Interfacial properties and our methods to compute them are described 
in section || The cluster algorithm for an efficient simulation of the Ising model 
with antiperiodic boundary conditions is explained in detail in section |j. Details of 
the data analysis and the numerical results are summarized in section [5]. Section | 
contains our conclusions. 

2 The model 

We consider a simple cubic lattice with extension L in x- and y-direction and with 
extension T = 2D + 1 in ^-direction. The lattice sites i = (i x ,i y ,i z ) have integer 
coordinates, and we adopt the convention that the z-coordinate runs from —D to 
+D. The Ising model is defined by the partition function 

Z= ]T eM-PH). (I) 
to=±i} 

The Hamiltonian if is a sum over nearest neighbor contributions, 

H = - k ij a i a j ■ ( 2 ) 

<i,j> 

The lattice becomes a torus by defining that the uppermost plane is regarded as the 
lower neighbor plane of the lowermost plane. An analog identification is done for 
the other two lattice directions. For the Ising spin field a we will use two different 
boundary conditions: Periodic boundary conditions are defined by letting kij = 1 
for all links < i,j > in the lattice. To define antiperiodic boundary conditions 
in ^-direction, we also set = 1 with the exception of the links connecting the 
uppermost plane (z = +D) with the lowermost plane (z = —D). These links carry 
an antiferromagnetic factor k^ = — 1. 

3 Interfacial properties 

In this section, we shall give a short account of important interfacial properties: 
surface width, surface tension and surface stiffness. 

3.1 Surface width 

For sufficiently large ft and large enough L, the imposure of antiperiodic boundary 
conditions forces the system to develop exactly one interface, a region where the 
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magnetization rapidly changes from a large negative value to a large positive value. 
(Situations where more than one interface can occur are discussed below.) An 
important property of an interface is its width. The definition of the interfacial 
width is not unique. We adopt the following definition: 

A magnetization profile for lattice planes perpendicular to the ^-direction is 
defined by 

M(i z ) = L- 2 J2 *i • (3) 

First we have to unwrap the torus on the real line. One then shifts the interface 
close to i z — by taking the following measure. An approximate interface position 
is defined as the value of i z where the absolute value of the magnetization profile 
takes its minimum. The whole Ising configuration is then shifted in z-direction such 
that this minimum comes close to i z = 0. Note that spins that pass the antiperiodic 
boundary at i z = ±D are to be flipped. 

We introduce an auxiliary coordinate z that assumes half-integer values (labeling 
positions between adjacent lattice layers perpendicular to the ^-direction), z takes 
values — D + 1/2, — D + 3/2, . . . ,D — 3/2, D — 1/2. A normalized magnetization 
gradient is introduced via 



P{z) 



M(z + \)-M{z- 1 -) 



2M b 

where is the positive bulk magnetization, i.e 



(4) 



M b = (L 2 T) 



-i 



E 



(5) 



For a given configuration of the spin field, the position of the interface is defined as 
the sum over zp(z). The interface width is defined as the expectation value 

w/ 2 = <Ep(^ 2 -(Ep(^) 2 >- (6) 

z z 

Especially on small lattices, fluctuations in the two bulk phases can deteriorate the 
results. Due to bubbles, p(z) can be accidentally large even far away from the 
interface position. Since such fluctuations contribute to the interface width with a 
weight proportional to the distance from the interface position, the true signal can 
disappear in the noise. One possibility to reduce noise that stems from fluctuations 
of bubbles in the bulk is to take the lattice extension T as small as possible. In the 
framework of a study based on a Metropolis algorithm, this approach is proposed 
in I|l3l EDI. However, one has to be careful not to disturb the free fluctuation of the 



interface (the properties of which we are interested in). In refs. |TI], T2\ it is proposed 
to implement a procedure to remove the bubbles before measurement of the actual 
magnetization profile. Note that when the bulk correlation length is small, then also 
the bubbles are small. One then can assume that the interface width changes little 
when one removes all bubbles. The procedure is as follows: All nearest neighbor 
pairs with a saturated bond are frozen together. (Note that a bond connecting top 
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and bottom layer of the lattice is saturated when the spins have different signs.) 
This defines a configuration of clusters. By flipping all spins in the largest cluster, 
all first order bubbles (bubbles which do not contain smaller bubbles) are completely 
removed. Iterating the procedure, one can quickly get rid of all bubbles from the 
configuration. We denote the interface width (measured as described above) on the 
bubble free configuration by Wq. 

3.2 Surface tension 

The surface tension a of a d-dimensional Ising model is defined by 



Here, Fi = — In Zi is the reduced free energy of the system with boundary condi- 
tions such that an interface perpendicular to the ^-direction is introduced at a fixed 
position. The boundary conditions of the system labeled by the subscript "0" are 
such that no interface is forced into the system. 

There are many possibilities to obtain estimates for the surface tension from 
Monte Carlo simulations on finite lattices, see e.g. refs. 0, 0, f|, [5], [7|, 

One has to do essentially with two sources of systematic errors that are distinct 
in nature but intimately related. Effects from too small L and effects from too small 
T. 

To minimize the finite size effects in z-direction we choose antiperiodic boundary 
conditions as described above. These boundary conditions do not fix the position of 
the interface: it can wander freely in z-direction and is less affected by the presence of 
a boundary compared to a system with fixed "H — " boundary conditions. However, 
it is still important that T is large compared to the width of the interface. 

The finite L-effect is as follows: For (3 C < (3 < j3 r , the interface is rough, which 
means that its thickness grows logarithmically with L. This means, that if we go to 
large interface areas we simultaneously have to increase T in order to avoid strong 
effects from confining a wildly fluctuating surface to a flat box. On the other hand, 
if we choose L too small, tunneling becomes more likely and there is a tendency 
that more than one interface will form. Note that in general there will be an odd 
number of interfaces for antiperiodic boundary conditions and an even number of 
interfaces in the case of periodic boundary conditions. For very small L the notion 
of an interface can even become meaningless. 

The finite size effects become the stronger the closer one approaches the bulk 
critical point where no interface survives the thermodynamic limit. This means that 
close to the bulk critical point one needs large and thick lattices to get systematic 
errors under control. 

Let us assume that there is only one interface in the system with periodic bound- 
ary conditions, and that there are no interfaces in the periodic system. Then the 
effect of the free motion of the interface in z-direction on the interfacial free energy 
amounts to add InT: 




(7) 



F s = F a . p . -F p . +lnT. 
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What happens to the surface free energy if there are several interfaces in the system? 
If one assumes that the interfaces do not interact which each other, one finds 

tanh(exp(-F s + InT)) = %^ . (9) 

If we resolve this equation with respect to F s we get 

F s = \nT-\n{hn{\±^&)). (10) 

^ 1 ~ ^a.p./^p. 

In general on has no direct access to the partition function in Monte Carlo simula- 
tions. (For not too large systems, the surface free energy can be obtained directly 
from a Monte Carlo simulation of a statistical ensemble that includes the boundary 
conditions as dynamical variables. These variables are updated using a modified 
cluster algorithm [|15|] .) 

In this paper we shall employ two methods to get estimates for the surface free 
energy. 

3.2.1 Surface free energy from integration over (3 

Note that the derivative of the free energy with respect to the coupling (3 is a well 
defined observable, 

dF 

Tr iH). (ID 
In the case of a single interface one therefore gets 

dF s 



dp 



H) a . p . - (H)p. . (12) 



Here, the expectation values are defined in the systems with periodic or antiperiodic 
boundary conditions, respectively. Let us introduce the abbreviation 

E s = {H) a .p.-{H)p.. (13) 

The surface free energy can then be obtained by integration over (3: 



F s ((3) = F s ((3 ) + f dl3' E S {13') . (14) 



Our approach is to compute by Monte Carlo simulation the surface energy for (3- 
values ranging from low temperatures around (3 = 0.6 up to the bulk critical region 
around (3 C . Note that we can integrate our data starting both from the hot and the 
cold side since the initial conditions for the integration are known in both cases: 
For large (3 we can employ a low temperature expansion for the interface tension by 
Weeks et al. (published in an article by Shaw and Fisher |l6|) to obtain the surface 
free energy of an interface at a fixed position. 

In the thermodynamic limit the surface tension vanishes in the high temperature 
phase, while it is finite in the low temperature phase. For finite systems the difference 
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F a . P . — F p strictly vanishes only at /3 — where only the entropy and not the energy 
enters the free energy. But F ap — F p will remain negligibly small until f3 comes 
close to (3 C . For our numerical purpose we set this point where the difference of 
the energies with periodic and antiperiodic boundary conditions exceeds a certain 
amount (which we define as the statistical error we can achieve in our simulations). 
Let us call this coupling (3\. 

Starting from f3\ where the difference F a/p . — F p . is negligible within the obtainable 
accuracy, we can integrate the energy differences to obtain the free energy for any 
f3. Note that in this case we have to take into account the existence of more than 
one interface for a certain range of the integration interval. In practice this means 
that we have to take InT as the integration constant: F s {(3\) = InT. 



3.2.2 Surface free energy from finite differences 

An alternative way to determine free energies is to add finite differences in the free 
energy from small intervals A/3. For sufficiently small A/3 and for not too large 
transverse lattice extension L we can get the change of the free energy directly from 
a single Monte Carlo simulation. It is easy to show that 

F{f3 + Af3) = F{f3) - ]n(exp(- APH)) P . (15) 

We put an extra subscript (3 here to make explicit that the expectation value is in 
the system simulated at inverse temperature (3. Results based on the usage of eq. 
( |15| ) can be easily checked for accuracy and consistency: by simulating at a certain 
point (3 one gets estimates for the free energies in a whole neighborhood of (3. Note 
that one can use negative A/5's as well. Now assume that we do another simulation 
at (3' > (3 not too far away from f3. Then we have two sets of estimates for the 
points between f3 and namely the ones from the simulation at (3 with positive 
A/3's and the ones from the simulation at (3' with negative A/3's. If all the results 
are consistent (within the statistical accuracy), we assume that the step from (3 to 
(3' was safe, and we proceed to the next larger /3-value. 



3.2.3 Surface tension from finite L data 

In this section we shall describe how we extract estimates for the surface tension a 
from the finite L data for the surface free energy F s . The fundamental definition of 
a as given in eq. (0) requires to actually perform the limit L — > oo. However, we 
observed that with very good precision the surface free energy behaves like 

F s = C s + a'L 2 (16) 

already for moderate surface extension L. It is therefore natural to identify the 
coefficient a' in eq. (|I6| ) with the surface tension a. 



Let us discuss how this definition of a surface tension on finite lattices relates to 
another one used in the literature. In ||, the surface tension is computed via the 
finite L behavior of the energy splitting due to tunneling E 0a , 

E 0a = Cexp(-aL 2 ). (17) 
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Note that the constant C is not identical with the constant C s introduced in eq. 
(16). However, there is an approximate relation between the two constants that can 



be derived by approximating the Ising system in a long cylinder by a 1-dimensional 
Ising model with a (3- value chosen according to 2(3 = F s . (This approximation 
assumes that the interfaces are sharply defined and do not interact with each other. 
Both conditions are fulfilled if (3 and L are large enough.) With the abbreviation 
v = exp(— 2(3) one finds 

E 0a = ln((l + v)/(l-v)). (18) 
For large (3 one has approximately E 0a « 2v, and thus 

2v = 2exp(-C, - ah 2 ) « Cexp(-crL 2 ) . (19) 

So we finally obtain the relation 

ln2-C s wlnC. (20) 

3.3 Surface stiffness 

In the theory of interfaces the surface stiffness coefficient k plays an important role. 
It is defined as follows. In generalization of the surface tension definition given in 
eq. (0) one defines 

a{9) = lim lim -, 1 1 — (Fj - F ) , (21) 

where by suitable boundary conditions in the system "J" a single interface is enforced 
that makes an angle 9 e.g. with the x-axis. Expanding for small inclination angle 9, 

a{9)/ cos(#) = (j(0) + a'(9) 9 + - « 9 2 + . . . , (22) 

one defines the stiffness coefficient k, 

K = (J (°) + ^l^o- ( 23 ) 

This coefficient plays an important role in the capillary wave model of rough inter- 
faces []17|. Roughly speaking, this model assumes that the surface dynamics of a 



rough interface is well described by a Gaussian model for surface "height" variables 
htXi, X 2 , Xd-i)- The model Hamiltonian is 

1 r d ~ 1 dh 

H cw = -J dX x dX 2 . . . dX d ^ K *(^) 2 > ( 24 ) 

~ * i=l 1 

where the K{ are the stiffness coefficients corresponding to inclinations of the interface 
with respect to the i'th lattice plane perpendicular to the cfth direction. In our case 
of a 3-dimensional Ising model on a simple cubic lattice, K\ — k% = k, and we define 

(3 cS = - ■ (25) 
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Capillary wave theory then says that the long distance properties should by encoded 
in a (d — l)-dimensional Gaussian model (massless free field theory) with partition 
function 




Long distance properties are most systematically studied via the block spin renor- 
malization group [|TI| . For the Gaussian model defined through eq. (p6|) one defines 
block spins $j as averages over cubic blocks I of size Ljf 1 : 

*i = L B (d - l) E hi ■ (27) 

Usually the renormalization group flow is described in terms of effective Hamiltoni- 
ans parameterized by effective coupling constants. For our purpose it is sufficient to 
consider expectation values of block spin observables which can be directly measured 
with the Monte Carlo method. We define the two quantities 

4° } = 4 E - ^) 2 ) > (28) 

' <I,J> 

where I and J are nearest neighbors in the block lattice, and 

Ai S = (4 E - M 2 ) , (29) 

1 [I,K] 

where / and K are next to nearest neighbors. I is the extension of the block lattice, 
i.e. I = L/Lb- For the Gaussian model, the A's can be computed exactly with the 
help of Fourier transformation. The results for a variety of lattice sizes are quoted in 
table 10. These values are computed for (3 = 1. Note that the results for arbitrary 
(3' can be obtained by just multiplying with (3' / (3. 

How do we now do the blocking for the Ising interface? Block spin "height 
variables" hi are defined as follows: Blocks / are defined as sets that are quadratic 
in x — y direction with extension Lb x Lb and that extend through the whole lattice 
in ^-direction. One block thus contains L\T lattice points. Now note that such a 
block can be regarded as if it were an Ising system in its own right. A magnetization 
profile and an interface position can be determined exactly as in the case of the full 
lattice. We define 

hi = interface position in block / . (30) 

Note that the blocked height variables can also be defined "with and without bub- 
bles" . Blocked observables for the Ising interface are introduced analogously to eqs. 
@ and (gg): 

4!r s) = 4 E fa - hi) 2 ) . (si) 
1 <i,j> 

where I and J are nearest neighbors in the block lattice, and 

4!r S) = <4Z (hi - h K ?) , (32) 

[I,K] 
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where / and K are next to nearest neighbor blocks. For a rough Ising interface, we 
define an effective coupling /3 eff as follows: 

^( Isin g) 

^ eff = r Um ^(0) ' ( 33 ) 

Of course, we expect that the so defined (3 e g does not depend on i or L In Monte 
Carlo studies we can not really do the infinite L^-limit. Instead, one has to convince 
oneself that the results have negligible finite size effects. This can, of course, only 
be checked within the given statistical accuracy. 



4 Cluster updating of Ising interfaces 



Straightforward application of the bulk cluster algorithm is not appropriate for 
a rough interface. The interface is correlated on all length scales, whereas the bulk 
correlation length is finite. Furthermore, competing interactions are induced by the 
interface, and cluster algorithms become inefficient for frustrated systems, see e.g. 



ref. . For the simulations reported on in this paper we used the interface cluster 
algorithm proposed in refs. |E| [H| and a slight modification of it, which is more 
suitable for simulations close to the bulk critical temperature. The motivation for 



these interface cluster algorithms stems from the VMR cluster |23| algorithm that 
allows efficient updating of 2-dimensional solid-on-solid (SOS) models. 

The basic idea of the algorithm is to map the system which is frustrated due to the 
antiperiodic boundary conditions stochastically onto a system without frustration. 
Then this auxiliary system is updated using a standard cluster algorithm [I9|, |20| . 
possibly with a small modification. 

The update of the auxiliary system fulfills detailed balance. For the entire update 
procedure, detailed balance can be implemented by making sure that the a priory 
probability of selecting a specific auxiliary system is not altered when the auxiliary 
system is updated. 

Let us now explain how we get the unfrustrated auxiliary Ising systems. The 
first step is to select a symmetry plane of the lattice perpendicular to the z-axis. 
This plane can be either put between two layers of sites or onto a layer of sites. Let 
us label the first alternative by "B" (for "between") and the second by "O" (for 
"onto"). 

(Strictly speaking, the torus is cut twice by the plane. For even T one has two 
B or two O positions, while for odd T there is one B and one O position. The 
algorithms discussed here will essentially make use of only one of these positions, 
depending on where the cluster growth is started.) 

We shall describe below how the reflection planes are selected stochastically with 
the right probabilities. 

Let us first consider the case that the plane is put between two adjacent layers. 
The plane divides the lattice into two parts A + and A_. Note that the points of 
A_ can be obtained from the points of A + by reflection with respect to the plane, 
denoted by i — > r(i). 
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The Ising variables S{ of the auxiliary system are now constructed to describe 
simultaneous flips of the spins <jj and a r (i)- For i 6 A + , one substitutes in the original 
Hamiltonian of the system 

Oi — > SiOi 

o"r(i) — > SjCr r (j) , (34) 

and determines the conditional probability distribution of the s-variables (the cr's 
kept fixed). For the partition function of the embedded Ising system one gets 

Z aux = eX P J ij S i S i > ( 35 ) 



where the couplings Jy in eq. (|35|) are given by 

Jij = P{kij<Ji<Jj + K(i)r(j)Vr{i)0r{j)) , (36) 

and the summation is over links in A + only. Note that the links connecting A + and 
A_ are not contained in the sum: Their energy is invariant under the possible spin 
flips. 

In the case of the reflection plane being identical with a layer of lattice sites 
one has to take into account that the sites on the plane are identical with their 
reflected partners. Hence the corresponding auxiliary spins control only one spin of 
the original system. This leads to a modification of the couplings J mn when both 
sites m and n are part of the reflection plane 

Jmn Pk mn <T m (T n , (37) 

while couplings connecting sites on the reflection plane with sites on the neighbor 
plane keep the form of eq. (pj|). 

In ref. |10| it is shown that there is no frustration in the auxiliary system. This 
means that the product of the J^ around an arbitrary closed path in A + is positive 
or zero. 

To make sure that the Boltzmann weight of the auxiliary configuration is equal 
to the Boltzmann weight of the original configuration we set the spins Sj = 1 before 
each update. After an update of the auxiliary system the and a r ^ are flipped if 
Si = —1. Otherwise they keep their old value. 

Let us now explain the choice of the reflection plane and the update of the 
auxiliary system. In the case of the modified (and simpler) version of the algorithm, 
one selects one of the possible reflection planes with uniform probability. Then one 
updates the auxiliary system by generating and flipping a single cluster using the 
delete probabilities 

f exp(- Jij(l + SiSj) if >0, . . 

Pdel > ij \ exp {+Jij{l - s iSj ) if J^ < . 1 J 

Since we want to update interface properties we start the cluster at a randomly 
chosen site on the reflection plane in the case of the O position and at a site next 
to the reflection plane in the case of the B position. 
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Let us label this version of the algorithm by "S" (for "simple"). This version is 
suitable when the width of the interface is of the order of the thickness of the lattice 
T. 

In the case of a well defined interface with a width small compared to T a 
more sophisticated procedure for the choice of the reflection plane is necessary. The 
corresponding version of the algorithm will be labeled with a "C" (for "cut"). As in 
the case of the VMR algorithm for SOS models, the crucial point is that the plane 
cuts the interface, thus dividing it into valleys and mountains which are candidate 
objects to be flipped. 

Again we shall describe the procedure first for reflection planes located between 
two adjacent layers of the lattice. For a given configuration of the Ising spins a, we 
consider the set of all unsaturated bonds that point in ^-direction. Let us denote 
this set by B. If there is a well defined interface, then obviously a large portion of 
B will be contained in the interface. The following prescription will therefore lead 
to a reasonable frequency of reflection planes that cut the interface: With uniform 
probability select a bond from the set B and take as reflection plane the plane that 
cuts this bond. 

We have to demonstrate that this choice respects the detailed balance condition: 
The first observation is that the selected bond will stay in B since its value cannot 
be changed by cluster flips. However, we have to make sure that this same bond is 
selected with exactly the same probability when the next cluster building is prepared. 
This can only be guaranteed if the size of B, i.e. the number of unsaturated links 
stays unchanged during the update. This can be achieved by using the following 
modified delete probabilities for bonds in ^-direction (for the other bonds the usual 
Swendsen-Wang choice defined in eq. ( |5%| ) stays in power) 

J 1 if J tj = , , . 

Pdel,ij j if J i:j + . 1 ' 

In order to properly select a reflection plane coinciding with one of the lattice 
layers, one proceeds as described above but afterwards shifts the reflection plane by 
±1/2, selecting each of the alternatives with probability 1/2, respectively. 

As in the case of the S algorithm we perform a single cluster update. We also 
start the cluster at a randomly choosen site on the reflection plane in the case of the 
O position and at a site next to the reflection plane in the case of the B position. 

To restore ergodicity, which is not satisfied by the deleting probabilities defined 
in eq. fl39|), we perform one Metropolis sweep after each update. 

To summarize, the Ising cluster algorithm comes in four brands that can be 
labeled by BS,OS,BC,OC. As was the case for the VMR algorithm for the SOS 
models, a mixture of two algorithms of type B and O is necessary to defeat critical 
slowing down. The S version performs better close to bulk critical point, and the C 
version is to be preferred if the interface is well defined. 

Let us close this section with the remark that one can show that the algorithm 
described in this section becomes identical with the VMR algorithm in the limit of 
infinite couplings in ^-direction, i.e. in the SOS limit. The BC version becomes the 
I-algorithm while the OC version corresponds to the H-algorithm of ref. [23]. 
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5 Data analysis and Monte Carlo results 



In this section we discuss in some detail how we evaluated our Monte Carlo data 
and summarize our results. 



5.1 Surface free energies and surface tension 



As described in section |3.2.1| , one of our methods to access the surface free energy 
is to integrate the surface energy over (3. Here we now describe how we did this in 
practice. 

We did simulations with antiperiodic boundary conditions in ^-direction on lat- 
tices with L — 8, 16, 32, 64 for /3-values ranging from the bulk critical region up to 
j3 = 0.6 which is deep in the low temperature domain. For many /5-values we made 
runs with different D to control the effects of a finite thickness of the lattice in 
^-direction. In total, we made more than 250 different simulations with antiperiodic 
boundary conditions. Typically, we made 10000 measurements of several quantities, 
separated always by 8 cluster updates (alternating of type O and B) and a single 
Metropolis sweep. A major part of the simulations was done on RISC workstations. 

A remark concerning the performance of the algorithms: The algorithm version 
C performed quite well for (3 > 0.25. The typical autocorrelation times (in units of 
the sweeps defined above) were of order one throughout. The type S worked well 
for all P, as long as the ratio of T and the surface width was sufficiently small. 

The simulations supplied us with a sufficiently dense grid of /^-values for the 
energies E a p . 

For most of the /3-values, we fortunately did not have to do extra simulations 
to access the energies with periodic boundary E p .. Instead we used the diagonal 
Pade approximation (order of nominator = order of denominator = 12) of the low 



temperature series by Bhanot et al. [25], cf. the appendix. By comparing with Monte 
Carlo simulations we found that this approximation is safe for a seizable range of 
(3- values, cf. table 1. In the table we quote the (3- values above which we used the 
Pade approximant throughout. For smaller (3 values we used the Cluster Monte 
Carlo method to determine E p ,. To give an impression of the data we display our 
results for the surface energies (divided by L 2 ) in figure 1. 

In order to do the integration over (3 we interpolated the data with the help 
of cubic splines which can easily be integrated over arbritrary intervals numerically. 
Estimates for the statistical error of the surface free energy were obtained as follows. 
For each of the /3-values, we have an energy value and an error bar. Note that the 
data for different (3 are statistically independent. We simulated a whole sequence 
of outcomes "energy as function of j3" by generating independent Gaussian random 
numbers centered around the Monte Carlo averages and with variances determined 
by the error bars, respectively. For each of these outcomes, a spline was generated 
and integrated. The error of the result of the integration (the free energy) was then 
obtained as the mean square deviation over this "data Monte Carlo" . 

We employed this method to do the integration over (3 starting from large (3 as 
well as from small (3. The integration from large (3 was always started at (3 = 0.6, 
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where the integration constant can be safely taken from the low temperature series 
for the surface tension. The "initial conditions" for the integration starting at small 
(3 were already described in section |3.2.1. 



The results for the free energies for L = 8, 16, 32 and 64 were then used to 
make fits with the ansatz eq. ( |T6| ) to obtain estimates for the surface tension a. 
Two typical fits are displayed in figure 2. We determined cr's both from the two 
different integration directions. A subset of our results is displayed in table 2. The 
X 2 's quoted in the last column show that the fits for the data from the integration 
started at small (3 have a significantly higher \ 2 - We also made fits with the L = 8 
data excluded, and they had a better x 2 /d.o.f. We conclude from this that close to 
the critical point the inclusion of larger lattices might be necessary to give reliable 
estimates for a. Our method would allow us to do this. 

With figure 3, we demonstrate how nicely the results from the two integration 
directions match within the error bars. The constants C s , however, show deviations 
indicating systematic effects which we do not have under control. This, of course, 
carries over to the approximate determination of the constant C defined in equation 



We also tried to estimate the critical exponent fi that governs the critical behavior 
of the surface tension. To this end we fitted our results for a according to the critical 
law a = <7oi M . Since our estimates for the cr's at different /5-values are strongly 
correlated we took the covariance matrix into account when doing the "data Monte 
Carlo" for the error estimates. We made two sorts of fits: Fits with the definition 
t = 1 — (3 C / (3 and fits with the definition t = (3/ (3 C — 1. We also varied the interval, 
over which the the (3 dependence of a was fitted. Our results are shown in tables 
3,4 and 5. The fits were always done using four different /3-values. Using more 
data points would not make very much sense since the data are correlated anyhow. 
Note however, that our statistical errors are nevertheless correct since our "data 
Monte Carlo" takes the covariances correctly into account. The comparison of the 
two different fits (using the two different definitions of t) clearly shows that there 
are systematic effects larger than the error bars: one still is not close enough to 
criticality. However we think that it is fair to say that our results are consistent 
with the value of \i ~ 1.26 expected from Widom scaling. The results for the critical 
amplitude <To show even stronger dependency on the type of the fit, and we can not 
say very much more that ln<7o is probably something between 0.2 and 0.4. 

The method to compute surface free energies by finite A/3-steps worked quite 
well. In tables 6 and 7 we present some of our results for the surface free energy on 
lattices with L = 8, 16,32 and 64. We there display the naive free energy F SjIl (L) 
obtained by assuming only a single interface, and the "improved" free energy F S ^(L) 
that is computed taking into account the presence of several interfaces. For small 
interface area and for j3 close to the critical point the difference between the two 
definitions becomes significant. We also determined estimates for the surface tension 
a from the surface free energies. The results are quoted in table 8. Again the fits 
have relatively large x 2 > an d in some cases discarding the L = 8 data changes the 
results beyond the statistical error. We again consider this as a warning that too 
small L's might lead to systematic errors in a. 
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5.1.1 The correlation length in the broken phase 

In order to study the behavior of the products £ 2 <r we tried also to extract the corre- 
lation length from the simulations of the systems with periodic boundary conditions. 
We define the magnetization of a time slice as 

St, = ^ E »t ■ (40) 

We define the following connected correlation function G(t) of the time slice mag- 
netizations: 

G(t) = (\S lz S tz+ t\) - hm <|^ +T |> . (41) 

For our practical purposes we used r = D which is the largest separation on the 
periodic lattices. For large t the correlation function G(t) is assumed to behave like 

G(t) oc exp(-|) ■ (42) 

(As a consequence of our definition of G(t), £ can not rigorously be considered as 
the inverse of a physical mass. For finite L the masses are split due to tunneling. 
However, the signal of tunneling is suppressed by the insertion of absolute values in 
eq. §!]).) 

We extracted estimates for £ from two subsequent values t and t + 1. One gets 

^(M + i)=- lnG(t) _; nG(t+1) . (43) 

The statistical error increases rapidly with increasing t. On the other hand the 
estimate for £ is spoiled by systematic errors if we choose a too small t. As a 
compromise we chose the t for our final estimate selfconsistently such that t > £ > 
t + 1. The results are quoted in table 9. One should remark that for all our data 
the lattice size was at least 12 times the correlation length. 



We also tried to extract the correlation length from low temperature series 26 



27]. Our analysis is based on the quantity 



A 2 = 6Xp( ~ rl) (44) 

where £ is the correlation length that controls the exponential decay of the correla- 
tion function. From the series for the correlation length given in |26| we derive the 
Taylor series for A 2 , 

A 2 = u 2 - u 3 + 10m 4 - 14m 5 + 93m 6 - 201u 7 

4731 o 33759 „ 115875 1n 

+ u s u 9 + u 

7 56 28 

295251 u 4847861 12 15341397 13 14 . 

- — u + ^r u — 286-* + °< tt >■ (45) 
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where u = exp(— 4/5). To get reliable numbers from this series one has to con- 
trol the singularities. We proceed as follows. We first multiply the series by 
(u — 0.41205) L25 (w + 0.336), where the first term stems from the critical singularity 
governed by the critical exponent v that we assumed to equal 1.25 here. Variation 
of the input for v in the range 0.625. ..0.63 led to changes in £ smaller than 3 per 
cent for (3 = 0.222. The second term stems from an unphysical singularity on the 
negative real line. Then we Taylor expand the result and perform a Pade aproxi- 
mation on it (type [7,6]). We then divide the result of the Pade approximation by 
[u — 0.41205) 5 / 4 (u + 0.336). Finally we determine £ by inverting eq. (|44]) . In figure 
4 we plot the quantity £t", with v = 0.625 as obtained from the measured £'s and 
of the series result. Within the statistical accuracy the results are nicely consistent. 



5.2 Surface stiffness 

In all simulations we measured the block spin correlation functions A^ sing ^ as defined 
in section [3l| for i = 1, 2 and / = 2, 4. We define auxiliary quantities 



.(Ising) 

tf£ = > ( 4 6) 



4? 



where the y4^°^'s are taken for L = 256 (cf. table 10), which is essentially the infinite 
block size limit. We thus get eight values that all (in the large L limit) should 
converge towards the same /3 e fj. In table 11 we show two examples for these eight 
values (for two different values of the lattice thickness T) at /3 = 0.24. The values 
for the different i, I are fairly consistent within the statistical accuracy. Closer to 
the critical point the estimates from the bubble free configurations are more stable. 
We therefore decided to use only the bubble free data for the determination of the 
/^-dependence of j3 e s. In tables 12 and 13 we present our estimates for (3 e ^. The 
values were determined by averaging over the two quantities (3^ with i — 1 and 
i = 2. The first two lines of table 12 shows that our results become unstable close 
to the critical point where /3 e ff diverges. For (3 = 0.43 and (3 = 0.45, both points are 
in the smooth phase, f3 e & decreases with increasing block size, in agreement with the 
Kosterlitz-Thouless picture of the roughening transition. This behavior is consistent 
with an infinite macroscopic stiffness for f3 > f3 r . 

In figure 5, we show our results for two combined quantities, namely £ 2 cr and £ 2 k. 
In the product £ 2 cx, the exponents fi and — 2u of the reduced temperature t should 
cancel, and we expect that this product should be fairly constant in a neighborhood 
of the critical point. The full line in figure 5 was obtained by combining our a's 
from the integration method with the correlation lengths as obtained from the Pade. 
Since we do not know the error of the Pade approximation of the low temperature 
series we base our error estimate for this quantity on our error bars for the measured 
correlation length as reported in table 9 and on the statistical errors on the surface 
tension a. We estimate the relative precision of our results for £ 2 <r to be around 
5 per cent for the smaller /3's, certainly better in the large f3 region. This takes 
into account statistical errors only. There might also be systematic errors (due to 
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too small L's) in the surface tension close to the critical point. They might be 
responsible for another 5 percent relative uncertainty. The points with error bars in 
figure 5 show the product £ 2 k. The plot shows that in the critical limit the surface 
stiffness becomes the same as the interface tension. This is a consequence of the 
restoration of rotational symmetry at the bulk critical point. In the theory of critical 
wetting, the following quantity plays an important role f28|: 

u{(3) = 1/(4tt< 2 ) . (47) 

In |H we find an estimate u{T/T c = 0.8) = 0.88. Our result u{T/T c = 0.8) = 
0.882(5) is in nice agreement with this prediction. In the limit (3 — > j3 r , Kosterlitz- 
Thouless theory states that f3 e g — > 2/n. We use the estimate for (3 T cited in the 
introduction, and find £(/3 e ff) — 0.3163 (from the Pade that here certainly is reliable). 
We then find a "KT- value" of £ 2 k which is 0.1572. 



5.3 Surface width 

In order to demonstrate the efficiency of the Hasenbusch-Meyer algorithm we redid 



the surface width computation of Mon et al. n] at (3 = (3 c /0.8 = 0.2771 on lattices 



of size L x L x 27, with L = 32, 64, ... , 512. Our results for W 2 and W$ are 
summarized in table 14. 

We performed fits of the data using the ansatz 

W 2 = const + ^ ln(L) (48) 

that is motivated by Kosterlitz-Thouless theory of a rough interface and that should 
become precise for large enough L (depending, of course, on (3). We performed 
several fits on subsets of the data. Our results for (3 e s are summarized in table 15. 



The errors were determined by a "data Monte Carlo" as described in section |5]T 
We conclude that the L = 32 data should not be included in the fit and estimate 
that (3 eS = 4.3(2). 

The result has to be compared with the f3 e g as obtained from the renormalization 
group quantities A introduced in section [373] . Because of the moderate statistics we 



used only quantities measured on the bubble free configurations. We computed the 
auxiliary quantities /3 1 ^, as defined in eq. (|46|), measured on bubble free configura- 
tions only. The quantities (3^ should converge to /3 e ff for large L. Our findings are 
summarized in table 16. One should not overemphasize the L = 512-results, which 
suffer a bit from poor statistics. Instead we focus on the L = 128 and on the L = 256 
results. Within the statistical accuracy these data are fairly compatible with the 
value for (3 eS obtained from the surface thickness fits. We interpret our results as 
a further confirmation that the long distance properties of the Ising interface are 
correctly described by a massless Gaussian dynamics. 
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6 Conclusions 



In this paper, we have presented a numerical study of the Ising interface in three 
dimensions over a wide range of temperatures. The method to obtain the surface 
free energies by integration over f3 requires many separate simulations but turns 
out to be practicable and useful. Inclusion of large interfaces is possible because of 
the usage of a cluster algorithm also for the simulations with antiperiodic boundary 
conditions. 

Our analysis of the surface tension indeed showed that closer to the critical point 
large interface extensions L are necessary to get reliable values (a has a tendency to 
come out too large when the lattices are too small). A high precision computation of 
o"o and \i seems difficult, the systematic effects from a too large reduced temperature 
t are strong. 

The large interfaces allowed us to study the infrared surface properties. We 
could confirm the massless Gaussian behavior in the rough phase and also extract 
the stiffness coefficient from the renormalization group behavior. 
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Appendix: Low temperature series for energy and 
surface tension 



For the convenience of the reader, we here cite the low temperature series for the 
energy density and for the surface tension. 

Bhanot et al. [^5| recently pushed forward the low temperature series for the 



energy to 24*^ order in the variable u = exp(— 4/3). One defines 

e = 3(1 — {aiCFj)) for i, j nearest neighbors. (49) 

The series for e is 

e = 12u 3 + 60m 5 -84m 6 + 420u 7 - 1056u 8 

+ 3756 u 9 - 11220u 10 + 37356 m 11 - 1 18164 u 12 + 389220 u 13 

- 1261932 m 14 + 4163592 m 15 - 13680288 m 16 + 45339000 u 17 

- 150244860 w 18 + 500333916 u 19 - 1668189060 u 20 
+ 5579763432 u 21 - 18692075820 u 22 

+ 62762602860 m 23 - 211062133044 m 24 + C(m 25 ) . (50) 

A low temperature series to order u 9 for the surface tension a was determined 
by Weeks et al. The coefficients can be found in a paper by Shaw and Fisher ]16| . 
The series is 



a = 2(3 -2u 2 - 2u 3 - 10u 4 -16w 5 

242 , 7 „ 4334 , 
- m 6 -150u 7 - 734u 8 u- 

O 9 



1 + 0(u 10 ) . 



(51) 
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Figure 1: The results for the surface energy per area as function of (3 for spatial 
lattice extension L = 8, 16, 32, 64. These are the data to be interpolated by splines 
and integrated over in order to determine the surface free energies. 



Figure 2: Two examples for our results for the surface free energy as function of L. 
These data are fitted with the law F s = C s + <jL 2 to determine the surface tension. 



Figure 3: With this figure, we demonstrate the consistency of the results for o 
obtained from the integration started at large (3 with that obtained from the in- 
tegration started at small (3. The left half of the figure contains only the "small 
j3" data, whereas the right part contains only the "large f3" data. Here they meet 
at (3 = 0.235. (The data from both methods cover the whole (3 range. We chose 
this specific presentation only to demonstrate the consistency.) For small (3, the a's 
obtained by integrating from below have smaller error bars. 



Figure 4: Comparison of our results for the correlation length in the broken phase 



with the Pade approximation of the low temperature series of ref. ^6|. The figure 
shows the result for the quantity £t", with v = 0.625. 



Figure 5: Plot of our results for the combined quantities £ 2 <r (full line) and £ 2 k 
(points with error bars). Details are explained in the text. 
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Table 1: Comparison of Monte Carlo results for energy per site (= 3(<jj<jj), with % 
and j nearest neighbors, periodic boundary conditions) with results from the Pade 
approximation of the low temperature series. In the last column we quote the (3 
value above which we consider the usage of the Pade approximation as safe. 
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Table 2: Results for the constant C s and for the surface tension a. Type 1 means: 
obtained by integration starting at large f3, type s means: obtained by integration 
starting at small (3. We also quote the logarithm of the constant C defined in 
section 3.2.3. The last column gives x 2 P er degree of freedom, averaged over the 
"data Monte Carlo" . 
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2.192(81) 
1.970(20) 


1.09 
7.45 


0.226 


1 

s 


2.770(83) 
2.539(22) 


0.01014(13) 
0.010318(44) 


2.077(83) 
1.846(22) 


1.33 
5.35 


0.227 


1 

s 


2.651(77) 
2.425(21) 


0.01310(12) 
0.013303(52) 


1.958(77) 
1.732(21) 


1.27 
4.65 


0.228 


1 

s 


2.554(76) 
2.324(21) 


0.01620(12) 
0.016406(55) 


1.861(76) 
1.631(21) 


1.39 
6.16 


0.229 


1 

s 


2.477(81) 
2.235(22) 


0.01938(11) 
0.019589(56) 


1.784(81) 
1.542(22) 


1.62 
7.43 


0.230 


1 

s 


2.388(73) 
2.151(21) 


0.02268(11) 
0.022863(68) 


1.695(73) 
1.458(21) 


1.84 
6.46 


0.231 


1 

s 


2.298(80) 
2.073(24) 


0.02606(11) 
0.026226(63) 


1.605(80) 
1.380(24) 


1.72 
5.82 


0.232 


1 

s 


2.232(76) 
1.999(26) 


0.02952(11) 
0.029684(82) 


1.539(76) 
1.306(26) 


1.27 
5.03 


0.233 


1 

s 


2.165(77) 
1.934(25) 


0.03305(11) 
0.033221(82) 


1.472(77) 
1.241(25) 


1.25 
5.27 


0.234 


1 

s 


2.097(73) 
1.877(26) 


0.03665(11) 
0.036808(89) 


1.404(73) 
1.184(26) 


1.20 
5.24 


0.235 


1 

s 


2.048(69) 
1.831(24) 


0.04025(10) 
0.040429(84) 


1.355(69) 
1.138(24) 


1.14 
5.25 


0.236 


1 

s 


1.999(72) 
1.785(31) 


0.04392(11) 
0.044077(95) 


1.306(72) 
1.092(31) 


1.13 
5.31 


0.237 


1 

s 


1.947(71) 
1.747(32) 


0.04760(10) 
0.047758(85) 


1.281(71) 
1.054(32) 


1.01 

5.68 


0.238 


1 

s 


1.909(75) 
1.707(35) 


0.05133(10) 
0.05149(10) 


1.216(75) 
1.014(35) 


1.15 
5.69 


0.239 


1 

s 


1.869(70) 
1.671(30) 


0.05509(10) 
0.05524(10) 


1.176(70) 
0.978(30) 


1.08 
5.35 


0.240 


1 

s 


1.840(66) 
1.636(33) 


0.05889(10) 
0.05905(11) 


1.147(66) 
0.943(33) 


1.11 

5.68 
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Table 3: Results for the critical exponent [i as obtained from fitting the logarithm of 
the surface tension with the critical law In a = In <t + fi In t. For this table we used 
the <r's from the integration starting at large (3. Type 1 is a fit with t — 1 — f3 c / (3, 
and type 2 is a fit with t = (3/ (3 C — 1. The fourth and fifth columns give the results 
using cr's obtained only from subsets of the data (L — 8 data excluded or L = 64 
data excluded, respectively). The fits were done using four equidistant j3- values. 
The numbers in square brackets give the average x 2 i n the "data Monte Carlo" , not 
divided by d.o.f. 



type of fit 


fit interval 


all lattice sizes 


-(L = 8) 


— (L = 64) 


1 
2 


0.224 - 0.2300 


1.274(16) [0.27] 
1.243(15) [0.57] 


1.276(19) [0.52] 
1.244(18) [0.90] 


1.266(31) [0.13] 
1.234(30) [0.08] 


1 
2 


0.224 - 0.2313 


1.270(12) [0.21] 
1.234(11) [0.48] 


1.275(17) [0.35] 
1.239(17) [0.72] 


1.257(28) [0.19] 
1.222(27) [0.20] 


1 
2 


0.224 - 0.2325 


1.278(9) [0.38] 
1.238(9) [1.18] 


1.282(12) [0.67] 
1.241(11) [1.38] 


1.252(22) [0.25] 
1.213(21) [0.33] 


1 
2 


0.224 - 0.2338 


1.277(8) [0.25] 
1.232(7) [1.02] 


1.281(11) [0.52] 
1.236(19) [1.33] 


1.259(21) [0.28] 
1.216(20) [0.34] 


1 
2 


0.224 - 0.2350 


1.276(8) [0.57] 
1.226(8) [0.96] 


1.281(11) [0.59] 
1.232(19) [1.24] 


1.248(18) [0.43] 
1.200(17) [0.25] 



Table 4: Same as table 3, however based on results for a obtained from the integra- 
tion starting at small f3. 



type of fit 


fit interval 


all lattice sizes 


-(L = 8) 


— (L = 64) 


1 
2 


0.224 - 0.2300 


1.252(5) [0.24] 
1.223(5) [1.07] 


1.268(7) [3.63] 
1.238(7) [7.41] 


1.240(7) [3.27] 
1.212(7) [1.21] 


1 
2 


0.224 - 0.2313 


1.253(4) [0.56] 
1.221(4) [1.30] 


1.266(5) [3.3] 
1.233(5) [8.5] 


1.242(7) [3.40] 
1.211(7) [1.31] 


1 
2 


0.224 - 0.2325 


1.254(4) [2.18] 
1.219(4) [0.87] 


1.268(5) [2.79] 
1.231(4) [8.58] 


1.244(7) [4.57] 
1.210(6) [1.72] 


1 
2 


0.224 - 0.2338 


1.257(4) [5.01] 
1.218(4) [0.86] 


1.270(5) [2.96] 
1.230(5) [9.34] 


1.246(6) [5.87] 
1.208(5) [1.40] 


1 
2 


0.224 - 0.2350 


1.261(3) [3.90] 
1.218(3) [1.83] 


1.269(3) [2.48] 
1.226(3) [11.11] 


1.249(5) [4.56] 
1.208(5) [0.70] 
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Table 5: Results for lno"o as obtained from fitting the logarithm of the surface 
tension In a with the critical law In a = lna + //mi. For this table we used the 
cr's from the integration starting at small j3. Type 1 is a fit with t = 1 — f3 c //3, and 
type 2 is a fit with t = f3/f3 c — 1. The fourth and fifth columns give the results 
using cr's obtained only from subsets of the data (L — 8 data excluded or L = 64 
data excluded, respectively). The fits were done using four equidistant /3-values. 
The numbers in square brackets give the average x 2 m the "data Monte Carlo" , not 
divided by d.o.f. 



type of fit 


fit interval 


all lattice sizes 


-(L = 8) 


— (L = 64) 


1 
2 


0.224 - 0.2300 


0.372(18) [0.24] 
0.232(17) [1.07] 


0.423(22) [3.63] 
0.281(22) [7.41] 


0.341(23) [3.27] 
0.206(22) [1.21] 


1 
2 


0.224 - 0.2313 


0.376(15) [0.56] 
0.223(15) [1.30] 


0.416(17) [3.30] 
0.261(17) [8.52] 


0.349(25) [3.39] 
0.239(24) [1.31] 


1 
2 


0.224 - 0.2325 


0.383(12) [2.17] 
0.216(11) [0.87] 


0.422(14) [2.79] 
0.253(13) [8.58] 


0.357(23) [4.57] 
0.198(21) [1.72] 


1 
2 


0.224 - 0.2338 


0.392(11) [5.01] 
0.214(11) [0.85] 


0.429(15) [2.96] 
0.248(14) [9.34] 


0.364(18) [5.87] 
0.192(16) [1.39] 


1 
2 


0.224 - 0.2350 


0.404(10) [3.90] 
0.213(10) [1.83] 


0.427(9) [2.48] 
0.234(9) [11.11] 


0.374(17) [4.56] 
0.189(16) [0.70] 



Table 6: Surface free energies obtained by the step-by-step method, for L = 8 and 
L = 16. We here display the naive free energy (subscript s) obtained by assuming 
only a single interface, and the "improved" free energy (subscript i) that is computed 
taking into account the presence of several interfaces. 



p 


F s , n (S) 


^,i(8) 


F s , n (16) 


^,i(16) 


0.222 


3.026(12) 


2.475(43) 


3.688(19) 


3.337(40) 


0.223 


3.060(14) 


2.589(41) 


3.882(24) 


3.684(36) 


0.224 


3.103(16) 


2.706(39) 


4.200(30) 


4.108(36) 


0.225 


3.156(19) 


2.827(39) 


4.649(32) 


4.614(40) 


0.226 


3.228(19) 


2.964(33) 


5.208(35) 


5.197(36) 


0.227 


3.311(20) 


3.103(30) 


5.849(37) 


5.846(37) 
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Table 7: Same as table 6, however for L = 32 and L = 64. 



p 


F S ,n(32) 


F s ,i(32) 


F S ,n(64) 


^,i(64) 


0.222 


4.340(28) 


4.272(32) 


5.82(11) 


5.82(11) 


0.223 


5.706(36) 


5.702(37) 


12.19(15) 


12.19(15) 


0.224 


7.841(43) 


7.841(43) 


21.90(17) 


21.90(17) 


0.225 


10.386(48) 


10.386(48) 


32.78(19) 


32.78(19) 


0.226 


13.179(52) 


13.179(52) 


44.21(20) 


44.21(20) 


0.227 


16.123(56) 


16.123(56) 







Table 8: Results of the fit of the form Fi m = C s + aL 2 for the data obtained by the 
step-by-step method. The last column gives x 2 P er degrees of freedom. 



p 


L's used 


C s 


a 


rVd.o.f. 


0.227 


8,16,32 


2.28(3) 


0.01356(7) 


10.7 


0.227 


16,32 


2.42(5) 


0.01338(9) 




0.226 


16,32,64 


2.61(4) 


0.01023(5) 


6.2 


0.225 


16,32,64 


2.76(4) 


0.00738(4) 


4.0 


0.224 


16,32,64 


2.96(4) 


0.00469(4) 


8. 


0.224 


32,64 


3.16(8) 


0.00457(6) 




0.223 


16,32,64 


3.20(4) 


0.00231(4) 


34. 


0.223 


32,64 


3.54(7) 


0.00211(5) 
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Table 9: Results for the correlation length £. t denotes the distance at which the 
correlation length is determined. 



p 


L 


T 


t 




0.2250 


64 


65 


3-4 


3.57(10) 


0.2275 


64 


33 


2-3 


2.45(7) 


0.2285 


32 


33 


2-3 


2.16(6) 


0.2300 


64 


33 


2-3 


2.04(6) 


0.2327 


32 


33 


1-2 


1.66(3) 


0.2350 


32 


33 


1-2 


1.48(3) 


0.2391 


16 


23 


1-2 


1.28(3) 


0.2400 


16 


23 


1-2 


1.22(3) 


0.2500 


16 


23 


0-1 


0.92(1) 


0.2600 


16 


11 


0-1 


0.76(1) 


0.2700 


16 


11 


0-1 


0.68(1) 


0.2800 


16 


11 


0-1 


0.61(1) 


0.2900 


16 


11 


0-1 


0.55(1) 


0.3000 


16 


11 


0-1 


0.52(1) 


0.3200 


8 


11 


0-1 


0.46(1) 


0.3400 


8 


11 


0-1 


0.40(1) 


0.3600 


8 


11 


0-1 


0.36(1) 



Table 10: Exact results for A\{ for (3 — 1 as obtained by Fourier transformation. L 
denotes the size of the 2-dimensional lattice. 



L 


,(o) 

AL,2 


,(o) 

^2,2 


,(o) 

Al,4 


Ao) 

^2,4 


4 





187500 





250000 






8 





136719 





187500 


0.293527 


0.380581 


12 





126721 





175926 


0.257225 


0.340481 


16 





123147 





171875 


0.243918 


0.326090 


24 





120565 





168981 


0.234147 


0.315663 


32 





119655 





167969 


0.230662 


0.311978 


48 





119002 





167245 


0.228148 


0.309333 


64 





118773 





166992 


0.227263 


0.308404 


96 





118609 





166811 


0.226628 


0.307739 


128 





118551 





166748 


0.226406 


0.307507 


256 





118496 





166687 


0.226191 


0.307282 
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Table 11: The auxiliary quantities for (3 = 0.24 on a 64 x 64 x T lattice, with 
(b) and without (nb) bubbles. 



% 


I 


(b) T = 19 


(nb) T = 19 


(b) T = 27 


(nb) T = 27 


1 


2 


16.88(44) 


16.98(47) 


16.94(47) 


16.59(46) 


2 


2 


16.82(52) 


16.85(52) 


16.27(56) 


16.06(55) 


1 


4 


16.82(20) 


17.18(20) 


17.64(22) 


16.86(21) 


2 


4 


16.79(23) 


17.04(23) 


17.35(25) 


16.79(24) 



Table 12: (3 eS for (3 = 0.24 . . . 0.32 as obtained from 2x2 and 4x4 blocking on 
lattices with L = 16, 32, 64. Only bubble free configuration were used. For each 

1 / 2 / 

value of /, the estimate for (3 eS was obtained by taking the average of (3 e ^ and (3 e g. 





L = 


= 16 


L = 


= 32 


L = 


= 64 


P 


I = 2 


1 = 4 


I = 2 


1 = 4 


I = 2 


/ = 4 


0.230 






62.3(36) 


73.6(89) 


45.5 (20) 


45.59(86) 


0.235 


28.84(61) 


28.30(43) 


28.67(42) 


28.65(20) 


24.50(72) 


25.16(23) 


0.240 


16.54(22) 


16.02(11) 


16.99(27) 


17.01(12) 


16.32(50) 


16.83(22) 


0.245 






12.24(19) 


12.42(8) 


12.46(32) 


12.47(13) 


0.250 


9.80(13) 


9.69(7) 


9.68(17) 


9.74(7) 


9.91(23) 


10.00(10) 


0.255 






8.26(14) 


8.12(5) 


8.02(17) 


8.11(7) 


0.260 


6.75(9) 


6.90(5) 


6.89(12) 


6.87(5) 


6.83(14) 


6.85(6) 


0.265 






5.78(10) 


5.95(5) 


5.88(11) 


5.92(5) 


0.270 


5.25(7) 


5.31(3) 


5.19(10) 


5.23(4) 


5.26(10) 


5.18(4) 


0.275 






4.57(8) 


4.63(3) 


4.65(8) 


4.65(4) 


0.280 






4.05(8) 


4.16(3) 


4.13(7) 


4.14(3) 


0.285 


3.81(6) 


3.89(2) 


3.83(7) 


3.81(3) 


3.79(7) 


3.80(3) 


0.290 






3.57(7) 


3.52(3) 


3.43(6) 


3.41(2) 


0.295 






3.21(5) 


3.19(2) 


3.17(5) 


3.19(2) 


0.300 


2.97(5) 


3.01(2) 


2.87(5) 


2.93(2) 


2.96(5) 


2.93(2) 


0.305 






2.71(4) 


2.73(2) 


2.72(4) 


2.70(2) 


0.310 






2.55(4) 


2.53(2) 


2.55(4) 


2.51(2) 


0.315 






2.37(3) 


2.36(1) 


2.33(4) 


2.34(2) 


0.320 






2.15(3) 


2.18(1) 


2.15(4) 


2.18(2) 
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Table 13: [3 cS for (3 = 0.325 . . . 0.45 as obtained from 2x2 and 4x4 blocking on 
lattices with L = 16, 32, 64. 





L = 


= 16 


L = 


= 32 


L = 


: 64 


p 


I = 2 


/ = 4 


I = 2 


1 = 4 


I = 2 


1 = 4 


0.325 


2.04(3) 


2.11(1) 


2.03(3) 


2.07(1) 


2.05(3) 


2.03(1) 


0.330 






1.93(3) 


1.95(1) 


1.94(3) 


1.92(1) 


0.335 






1.80(3) 


1.81(1) 


1.79(3) 


1.82(1) 


0.340 






1.71(2) 


1.73(1) 


1.72(3) 


1.70(1) 


0.345 






1.58(2) 


1.62(1) 


1.58(2) 


1.60(1) 


0.350 


1.57(3) 


1.60(1) 


1.53(2) 


1.53(1) 


1.51(2) 


1.52(1) 


0.355 










1.40(2) 


1.42(1) 


0.360 










1.32(2) 


1.35(1) 


0.365 










1.23(2) 


1.27(1) 


0.370 


1.20(1) 


1.28(1) 


1.19(1) 


1.22(1) 


1.19(2) 


1.20(1) 


0.375 










1.13(2) 


1.13(1) 


0.385 










1.00(1) 


1.01(1) 


0.390 


0.94(1) 


1.00(1) 


0.93(1) 


0.96(1) 


0.94(1) 


0.94(1) 


0.405 










0.74(1) 


0.74(0) 


0.410 


0.70(1) 


0.76(1) 


0.66(1) 


0.69(1) 






0.430 


0.44(1) 


0.50(0) 


0.35(1) 


0.39(0) 


0.22(0) 


0.27(0) 


0.450 


0.24(0) 


0.31(0) 


0.13(0) 


0.17(0) 


0.05(0) 


0.08(0) 



Table 14: Squared surface width W 2 (measured on configurations with bubbles) and 
Wq (measured after removal of the bubbles) at /3 — /3 c /0.8 on lattices L x L x 27. 
stat denotes the number of measurements. Between two subsequent measurements 
we always performed eight single cluster updates (alternating type O and B) and 
a single Metropolis sweep. For L = 256 and L = 512 we show results of two 
independent runs. 





L = 32 


L = 64 


L = 128 


L = 256 


L = 512 


w 2 


2.341(42) 


2.889(42) 


3.365(29) 


3.855(35) 
3.811(33) 


4.313(60) 
4.335(64) 


w 2 


3.094(9) 


3.619(16) 


4.098(20) 


4.573(35) 
4.542(30) 


5.059(66) 
5.069(58) 


stat 


10000 


3000 


2430 


550 / 820 


187/260 



32 



Table 15: Fit results for /3 e ff ■ The numbers in square brackets give the average \ 2 
per degree of freedom in the "data Monte Carlo" . 



data used 


all 


all - (L — 32) 


all - (L = 32, 64) 


all - (L = 512) 


from W 2 


4.45(12) [1.33] 


4.30(16) [1.05] 


4.30(23) [1.22] 


4.44(14) [1.87] 


from Wq 


4.48(6) [2.29] 


4.3(1) [1.31] 


4.3(2) [1.70] 


4.48(6) [3.62] 



Table 16: ff^ for (3 = /3 c /0.8 as obtained from 2x2 and 4x4 blocking on lattices 
with L = 32, . . . ,512. The quantity (3^ is defined in the text. For L = 256 and 
L = 512 we show results of two independent runs. 



L 


Peflf 


PcS 


PcS 


PcS 


32 


4.50(8) 


4.41(9) 


4.46(3) 


4.45(4) 


64 


4.7(2) 


4.74(20) 


4.50(6) 


4.48(7) 


128 


4.04(18) 


4.03(21) 


4.433(82) 


4.434(95) 


256 


4.26(37) 
3.9(3) 


4.07(39) 
3.8(3) 


4.37(17) 
4.30(15) 


4.35(17) 
4.29(16) 


512 


4.65(62) 
5.1(6) 


3.16(64) 
5.9(8) 


5.02(35) 
4.25(21) 


4.62(31) 
4.6(3) 
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